Write a custom CUDA kernel to optimize Poly1FocalLoss (based on PolyLoss ICLR 2022 paper).

Formula: Loss = FL + epsilon * pow(1 - Pt, gamma + 1)
Where FL (Focal Loss) = -alpha_t * pow(1 - Pt, gamma) * log(Pt).
Pt is the probability of the ground truth class: p if target=1, 1-p if target=0.
p is sigmoid(logit).

Problem Analysis:
1. Memory Bandwidth: The standard implementation involves a long chain of element-wise operations: sigmoid, comparisons, subtractions, powers, logs, multiplications. Each step allocates and reads/writes global memory.
2. Mathematical Redundancy: Calculating the Focal term and the Poly-1 term separately involves computing the power function twice.

Optimization Strategy: Fully Fused Element-wise Kernel with Vectorization

1. One-Thread-per-Element: The loss is strictly element-wise. We launch a grid to handle N * C elements.

2. Vectorized Loads (float4): Load 4 logits and 4 targets (floats) at once into registers. This reduces memory instructions by 75%.

3. In-Register Math Simplification:
   Factorize the common term (1 - Pt)^gamma:
   Loss = (1 - Pt)^gamma * ( -alpha_t * log(Pt) + epsilon * (1 - Pt) )
   This reduces two expensive `pow` calls to a single `pow` call and a simple multiply-add.

4. Fast Math Intrinsics: Use `__expf`, `__logf`, `__powf` for maximum throughput on GPU special function units.

5. Reduction: The kernel outputs element-wise loss. The final reduction (mean/sum) is handled by ATen primitives in the C++ wrapper.
  
Here's an example to show you the syntax of inline embedding custom CUDA operators in torch: The example given architecture is:   
  
```python
import torch
import torch.nn as nn
import torch.nn.functional as F

BATCH_SIZE = 128
NUM_ANCHORS = 64
NUM_CLASSES = 64
SHAPE = (BATCH_SIZE, NUM_ANCHORS * NUM_CLASSES)

GAMMA = 2.0
ALPHA = 0.25
EPSILON = 1.0

class Poly1FocalLoss(nn.Module):
    """
    Optimized Poly1FocalLoss Implementation.
    L = L_fl + epsilon * (1-Pt)^(gamma+1)
    """
    def __init__(self, gamma=2.0, alpha=0.25, epsilon=1.0, reduction='mean'):
        super(Poly1FocalLoss, self).__init__()
        self.gamma = gamma
        self.alpha = alpha
        self.epsilon = epsilon
        self.reduction = reduction

    def forward(self, logits: torch.Tensor, targets: torch.Tensor) -> torch.Tensor:
        # 1. 计算 BCE Loss (即 -log(Pt))
        ce_loss = F.binary_cross_entropy_with_logits(logits, targets, reduction='none')

        # 2. 计算预测概率 p
        p = torch.sigmoid(logits)

        # 3. 计算 (1 - Pt)
        p_t_complement = (p - targets).abs()

        # 4. 计算 Alpha_t
        alpha_t = torch.where(targets == 1, self.alpha, 1.0 - self.alpha)

        # 5. 计算 Poly1 Focal Loss
        modulating_factor = p_t_complement.pow(self.gamma)
        
        loss = modulating_factor * (alpha_t * ce_loss + self.epsilon * p_t_complement)

        if self.reduction == 'mean':
            return loss.mean()
        elif self.reduction == 'sum':
            return loss.sum()
        return loss

class Model(nn.Module):
    def __init__(self, gamma=2.0, alpha=0.25, epsilon=1.0, reduction='none'):
        super(Model, self).__init__()
        self.loss_fn = Poly1FocalLoss(gamma=gamma, alpha=alpha, epsilon=epsilon, reduction=reduction)
    
    def forward(self, logits, targets):
        return self.loss_fn(logits, targets)

def get_inputs():
    logits = torch.randn(SHAPE, dtype=torch.float32)
    targets = torch.randint(0, 2, SHAPE, dtype=torch.float32)
    
    return [logits.contiguous(), targets.contiguous()]

def get_init_inputs():
    return [GAMMA, ALPHA, EPSILON, 'none']